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Introduction 


On-chip  capillary  electrophoresis  was  first  developed  in  the  early  1990’s  [1]  and  remains  an 
important  separation  technique  for  microfluidic  systems  today.  The  motivation  to  shrink  separation 
systems  comes  from  their  increased  portability,  reduced  consumption  of  expensive  reagents,  and  more 
feasible  parallel  analysis  of  separations.  Since  the  separation  efficiency  increases  with  the  separation 
length,  serpentine  channels  have  become  a  preferred  design  topology  because  of  their  high  channel 
density.  However,  introducing  turns  into  the  separation  channel  creates  additional  sources  of  variance 
from  the  new  electric  field  structure  and  geometric  path  length  differences. 

These  additional  sources  of  dispersion  have  been  previously  modeled  in  [2, 3, 4, 5].  The  problem 
with  each  of  these  models  is  that  they  are  based  on  the  assumption  of  a  linear  electric  field  structure  in  a 
single-turn,  which  excludes  many  higher-order  effects.  This  work  shows  a  new  model  of  the  high  Peclet 
variance  for  entire  serpentine  channels  using  the  full  structure  of  the  electric  fields  described  by  complete 
solutions  to  Laplace’s  equation.  This  model  is  based  on  closed-form  algebraic  expressions  for  rapid 
computation  and  possible  use  in  gradient-based  optimization.  It  is  a  functional  model  type,  as  described 
in  [6],  for  serpentine  channels  of  an  arbitrary  number  of  turns.  The  new  variance  description  can  be  used 
in  a  component-based  model  for  lab-on-a-chip  systems  that  consist  of  an  injector,  serpentine  separation 
channel,  and  detector,  which  is  a  common  separation  system  topology  [7]. 
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1.  Background 


Electrophoresis  is  the  method  used  for  species  transport  in  separation  microchannels  for  this 
work.  There  are  four  main  types  of  electrokinetic  phenomena  [8,9,10]:  Electrophoresis,  electroosmosis, 
sedimentation  potential,  and  streaming  potential.  The  first  two  are  the  main  types  of  transport  used  in 
electrokinetic  separation  systems.  This  section  will  briefly  introduce  double  layer  theory,  and  then 
discuss  its  application  to  electroosmosis  and  electrophoresis  in  order  to  form  the  basic  set  of  equations 
governing  the  dynamics  of  the  separation  systems  in  this  report. 

1.1  The  Double  Laver 

In  general,  when  two  phases  form  an  interface  (such  as  the  liquid  buffer  and  solid  insulating 
wall),  a  spontaneous  separation  of  charge  occurs.  The  separation  of  the  charge  into  two  components  at 
the  dual-phase  interface  is  what  gives  the  double  layer  its  name.  As  stated  by  Hunter  [10],  there  are 
several  reasons  for  such  a  separation: 

1 .  differences  between  the  two  phases’  affinity  for  electrons, 

2.  differences  between  the  two  phases’  affinity  for  ions, 

3.  ionization  of  surface  groups,  and 

4.  entrapment  of  fixed  charge  in  one  phase. 

The  specific  details  of  each  of  the  four  processes  above  are  quite  involved  and  Hunter  describes 
the  specifics  in  [10].  These  four  processes  are  responsible  for  the  creation  of  a  surface  charge  on  the 
channel  walls  and  the  formation  of  the  double  layer  between  the  wall  and  the  buffer.  This  surface  charge 
causes  dissociated  ions  or  other  charge  in  the  buffer  of  opposite  polarity  to  accumulate  in  a  double  layer 
near  the  wall,  as  seen  in  figure  1 .  In  general,  there  are  ions  of  both  types  of  charge  in  an  electrolyte  buffer 
solution.  However,  the  concentration  of  counterions  (ions  shielding  the  surface  charge)  is  generally  much 
greater  than  that  of  co-ions  (ions  of  charge  similar  to  the  wall’s  surface).  The  positive  and  negative 
charges  in  the  mobile  half  of  the  double  layer  (the  diffuse  portion  that  exists  within  the  liquid  buffer)  are 
described  by  Boltzmann’s  equation  [10]: 
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where  n±  is  the  local  positive  or  negative  ion  concentration,  n±°  is  the  ion  concentration  in  the  bulk,  zi  is 
the  valency  of  ion  /,  qe  is  the  charge  of  an  electron,  O  is  the  electric  potential  that  goes  to  zero  as  n±  goes 
to  n±° ,  k  is  Boltzmann’s  constant,  and  T  is  the  temperature.  When  <E>  is  positive  as  with  a  positively 
charged  surface,  n+  <  n+° ,  i.e.,  positive  ions  are  repelled  from  the  double  layer,  and  n.  >  n°  which  means 
negative  ions  are  attracted  to  the  double  layer. 

1.2  Electroosmosis 

For  electroosmosis,  when  an  electric  field  is  applied  axially  along  the  channel,  the  high 
concentration  of  counterions  in  the  double  layer  feels  a  force  and  begins  to  move.  This  concentration  of 
counterions  dominates  any  opposing  drag  of  co-ions.  The  net  effect  drags  all  subsequent  layers  adjacent 
to  the  double  layer  and  this  slip-boundary  condition  (which  is  a  valid  analogy  when  the  double  layer  is 
small  compared  to  the  channel)  quickly  forms  a  plug  flow  in  the  channel  [9].  Since  the  buffer  is  the 
system  component  that  moves,  this  type  of  transport  does  not  require  that  the  species  itself  be  charged. 


loose  counter-ions  in  the  ‘shield’  slide,  creating  a  plug  flow 
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Figure  1  -  (a)  Double  layer  formed  on  channel  walls  in  electroosmotic  separation,  (b)  Motion  of  the  double 
layer  and  resulting  plug  flow. 

The  double  layer  is  not  a  well-defined  boundary.  The  composing  ions  are  subject  to  thermal 
excitations  and  other  free  particle  interactions,  which  cause  them  to  diffuse  and  create  a  diffusion 
boundary  around  the  charged  surface.  One  way  to  characterize  the  double  layer  is  to  apply  Poison’s 

equation,  V20  —  pie,  within  the  layer  and  solve  for  the  potential.  There  will  be  a  certain  potential  a 
small  distance  away  from  the  charged  surface  where  the  particles  in  the  double  layer  will  be  able  to  flow 
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past  the  surface.  At  this  shear  surface,  the  potential  is  called  the  zeta  potential,  £,  [8,9,10].  At  a  certain 
distance  beyond  the  zeta  potential  point,  there  will  be  a  location  where  the  double  layer  potential  will 
decay  to  l/e  of  its  original  value.  Physically,  this  He  distance  is  the  point  where  the  potential  energy 
approximately  equals  the  thermal  energy  of  the  shielding  ions  [8].  This  distance  is  called  the  Debye 
length,  and  is  frequently  used  to  characterize  the  size  of  the  double  layer. 

1.3  Electrophoresis 

For  electrophoresis,  the  transport  method  used  in  this  work,  the  physics  are  very  similar  to  that  of 
electroosmosis.  When  charged  particles  are  placed  into  a  buffer,  they  begin  to  attract  charged  ions  from 
the  buffer  to  form  a  double  layer  as  in  figure  2. 


loose  counter-ions  in  the  ‘shield’  swim 


Figure  2  -  (a)  Double  layer  around  charged  particles  in  electrophoretic  separation,  (b)  Motion  of  double 
layer  and  resulting  net  particle  motion. 


When  an  electric  potential  is  applied  to  the  channel,  the  double  layer  around  these  charged 
particles  feels  an  electric  force.  Since  the  normal  component  of  the  electric  field  within  the  double  layer 
is  generally  much  stronger  than  that  of  the  applied  electric  field,  only  the  tangential  component  of  the 
applied  electric  field  has  an  effect  on  the  particle  [9].  This  causes  the  mobile  charge  in  the  double  layer  to 
“swim”,  which  exerts  a  viscous  force  on  the  buffer,  and  this  propels  the  charged  particle.  For  example,  in 
figure  2,  a  negative(positive)  potential  on  the  left  and  a  positive(negative)  potential  on  the  right  would 
cause  the  shielding  positive  charge  in  the  double  layer  to  migrate  to  the  left(right),  which  would  create  the 
viscous  forces  with  the  buffer  to  cause  the  charged  particle  to  “swim”  to  the  right(left).  From  this  point  of 
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view,  it  becomes  apparent  that  electrophoretic  and  electroosmotic  flows  are  complements  of  each  other, 
and  their  analysis  is  mostly  interchangeable  [8]. 

1.4  Equations  for  the  Electrophoretic  Dynamics  of  Separation  Systems 

When  looking  at  electrophoretic  motion  of  more  than  just  a  single  particle  and  assuming  that  the 
double  layers  have  a  small  Debye  length  and  remain  undistorted  in  the  flow,  and  that  the  flow  is  inertia 
free,  the  equations  describing  the  dynamics  of  the  species  flow  reduce  to  three.  First,  the  electric  field  in 
the  channel  is  described  by  Laplace’s  equation: 

V2<t>  =  0,  (2) 

since  for  a  small  Debye  length,  the  distance  away  from  the  charged  particle  before  it  appears  completely 
shielded  and  charge  neutral  is  not  far.  So  for  any  Debye  length  much  less  than  the  particle  radius,  the 
particle  appears  uncharged  at  any  length  scale  characteristic  of  the  channel  size.  Second,  the  inertia-free 
motion  of  the  species  is  related  to  the  electric  field  through  its  mobility,  which  is  described  by  the 
Helmholtz-Schmolochowski  equation: 

u  =  =  /. lE ,  (3) 

Mv 

where  £  is  the  permittivity  of  the  buffer,  £  is  the  zeta  potential  of  the  particles,  jUv  is  the  viscosity  of  the 
buffer,  and  E  is  the  electric  field.  Note  that  the  combination  of  £,  £  and  jUv  in  (3)  describes  the  mobility  of 
the  particles.  Finally,  with  the  inertia-free  effects  of  diffusion,  the  field-coupled  motion  of  the  species 
concentration  is  described  by  the  advection-diffusion  equation: 

—  +  m-Vc  =  DV2c  (4) 

dt 

The  remainder  of  this  thesis  will  use  the  previous  three  equations,  (2),  (3),  and  (4),  to  analyze  and 
produce  models  for  serpentine  separation  channels.  In  the  next  chapter  of  this  report,  the  field  structure 
will  be  derived  from  (2)  and  analyzed  for  use  in  (3).  Then,  chapter  3  will  explore  the  effects  of  the 
velocity  field  from  (3)  on  the  variance  of  the  concentration  distributions  of  analyte  bands  described  by 
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(4).  All  of  the  models  introduced  in  the  following  chapters  of  this  thesis  are  verified  within  their  chapter 
of  origin.  Finally,  chapter  4  will  provide  a  conclusion  for  this  work. 

Within  the  scope  of  serpentine  channel  electrokinetic  separation  systems,  there  are  several  types 
of  designs  possible.  The  use  of  modified  turn  geometries  is  proposed  in  [2,3],  however,  this  thesis  will 
focus  on  uniform  180°  turn  geometries  since  they  are  still  very  common. 
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2. 


Electric  Field  Structure  in  Uniform  Curved  Microchannels 


This  chapter  will  introduce  models  of  the  electric  field  for  serpentine  channels,  which  gradually 
relax  the  assumptions  made  in  [2, 3, 4, 5].  These  models  will  be  verified  with  finite  element  software 
where  necessary  as  they  are  developed.  The  next  chapter  will  use  these  new  models  to  create  more 
accurate  models  of  the  band  variance  in  the  separation  channels. 

Previous  work  [2, 3, 4, 5]  makes  two  assumptions  about  the  electric  field  in  their  single  turn 
dispersion  models.  The  first  assumption  is  that  the  electric  field’s  magnitude  changes  linearly  in  the  radial 
direction  within  the  turn,  which  is  an  approximation  of  a  field  whose  magnitude  varies  as  1/r  and  is 
circumferentially  directed.  The  second  assumption  is  that  the  field  changes  abruptly  from  a  constant 
value  in  the  straight  section  to  the  linear  circumferentially  directed  field  in  the  turn.  These  two 
assumptions  result  in  a  simple  linear  skew  of  the  band  after  the  turn,  whose  variance  is: 


where  skew=26w ,  which  represents  the  length  shown  in  figure  3,  the  turn  angle,  0. \  equals  ;rfor  180°  turns, 


and  w  is  the  width  of  the  channel,  [11].  An  interpretive  comparison  of  the  resulting  band  shape  with  and 
without  the  two  assumptions  made  in  [2, 3, 4, 5]  is  shown  in  figure  3. 


Figure  3  -  Actual  1/r  velocity  profile  shown  on  the  left,  linear 
approximation  shown  on  the  right.  As  the  width  of  the  channel  increases, 
the  accuracy  of  the  linear  approximation  decreases. 

The  assumptions  of  a  linear  circumferentially  directed  electric  field  (and  the  resulting  velocity 
profile)  with  the  abrupt  transition  in  the  electric  field  between  the  straight  leg  and  curved  turns  make  the 
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analysis  of  microchannel  dispersion  much  easier.  However,  such  analysis  is  only  valid  for  channels  with 
a  center  radius,  rc,  much  greater  than  the  width  of  the  channel.  A  picture  of  the  electric  field  components 
in  a  channel  with  a  w/rc  of  1.0,  along  axial  lines  at  w! 4,  wl 2,  and  3w/4,  is  seen  in  figure  4. 


Abrupt  (linear):  Plot  of  Component  Tangent  to  Path  Abrupt  (linear):  Plot  of  Component  ’Normal  to  Path 
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Figure  4  -  Plot  of  the  axial  and  transverse  components  of  the  abrupt-transition  linear  electric  field.  Note, 
there  are  no  transverse  components  as  shown  in  the  figure  on  the  right  for  the  abrupt  description.  This  is 
for  a  channel  with  a  width-to -center  radius  ratio,  w/rc ,  of  1.0. 


2.1  Abrupt-Transition  1/r  (1  -Assumption)  Field  Model 

This  section  will  reduce  the  number  of  assumptions  made  in  [2, 3, 4, 5]  by  not  linearizing  the 
electric  field  in  the  curve.  Since  the  abrupt-transition  assumption  is  still  made,  it  is  apparent  that  this 
description  will  also  be  missing  information  about  the  true  structure  of  the  electric  field,  as  seen  in  the 
equipotential  lines  of  figure  5a  and  5b.  However,  since  the  abrupt-transition  1/r  field  structure  contains 
more  information  than  the  simple  linear  models  and  is  much  easier  to  manipulate  than  a  full  Laplacian 
solution,  its  derivation  is  provided  in  this  section.  Note,  that  further  on  in  this  report  the  abrupt-transition 
1/r  field  structure  is  commonly  referred  to  only  as  an  ‘abrupt’  or  ‘abrupt-transition’  field  structure,  since 
the  abrupt  transition  is  its  only  assumption.  This  is  different  from  the  abrupt-transition  linearized  field 
structure. 
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Figure  5  -  (a)  Abrupt-transition  1/r  field  equipotential  lines  (b)  Laplacian  electrostatic  transition 
field  equipotential  lines,  (c)  Electric  field  lines  corresponding  to  (b) 


By  taking  the  line  integral  of  the  electric  field  along  the  channel  axis,  setting  it  equal  to  the 
voltage  applied  between  the  two  leg  ends,  and  requiring  the  conservation  of  current  between  the  straight 
and  curved  sections,  the  following  field  structure  is  derived: 


,  VJnjbUQ 
y  2l\n(b/a)  +  w/r 


(6a) 


wVo 

r(2/ln(Z?  /  a)  +  w7l)  ’ 


(6b) 


where  w,  r,  Z,  b ,  a ,  and  V0  are  the  dimensions  and  applied  voltage  as  seen  in  figure  5b.  These  results  show 
a  constant  field  magnitude  in  the  axial  direction  within  the  straight  legs,  as  seen  in  (6a),  and  a  1/r 
magnitude  field  in  the  axial  direction  in  the  curved  sections,  as  seen  in  (6b). 

Figure  6  visualizes  the  abrupt-transition  1/r  field  components  as  described  by  (6a)  and  (6b)  along 
axial  lines  at  w/4,  wl 2  and  3w/4.  Note  first  that  this  abrupt  description  has  no  transverse  components, 
which  completely  ignores  the  bending  of  the  equipotential  lines  around  the  intersection  of  the  straight  and 
curved  sections  seen  in  figure  5b.  Second,  note  the  nonphysical  transition  from  the  constant  magnitude 
field  in  the  straight  legs  to  the  1/r  circumferential  field  in  the  curved  section.  These  missing  components 
affect  the  variance  of  an  analyte  band  traveling  under  the  influence  of  the  electric  field  through  a  channel, 
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as  seen  in  the  next  chapter,  and  will  be  accounted  for  in  the  next  section  with  a  full  Laplacian  description 
of  the  electric  field  structure. 


Abrupt  (1/r):  Plot  of  Component  Tangent  to  Path 
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Figure  6  -  Plot  of  the  axial  and  transverse  components  of  abrupt-transition  1/r  electric  field.  Note,  there  are  no 
transverse  components  as  shown  in  the  figure  on  the  right  for  the  abrupt-transition  1/r  description  either.  This  is  for  a 
channel  with  a  width  to  center  radius  ratio  of  1.0. 


2.2  Continuous-Transition  Laplacian  (0-Assumption)  Field  Model 

The  full  Laplacian  solution,  as  seen  in  the  equipotential  lines  of  figure  5b,  has  more  structure  than 
the  abrupt-transition  1/r  solution  due  to  the  curvature  of  the  equipotential  lines  before  and  after  the 
straight  section’s  intersection  with  the  curve.  Additional  solutions  from  Laplace’s  equation  superposed 
with  the  abrupt-transition  solutions  from  the  previous  section,  which  are  valid  solutions  to  Laplace’s 
equation,  provide  this  additional  structure,  and  a  complete  model  that  agrees  well  with  full  numeric  PDE 
solutions  of  the  field  structure. 

To  derive  the  additional  solutions  to  Laplace’s  equation,  symmetry  is  utilized  and  only  half  of  the 
channel  is  analyzed,  as  seen  in  figure  7.  The  analysis  on  this  half  of  the  channel  requires  splicing  of  the 
solutions  in  the  straight  and  curved  sections,  which  can  be  achieved  with  suitable  boundary  conditions. 
To  obtain  the  correction  fields,  only  the  solutions  to  Laplace’s  equation  that  decay  away  from  the 
intersection  line  will  be  sought,  since  the  transition  that  is  apparent  in  the  equipotential  lines  of  figure  5b 
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fades  away  leaving  the  abrupt-transition  solutions  in  the  curve  and  straight  sections  far  from  the 
intersection  line,  whose  solutions  we  already  have  in  (6a)  and  (6b). 


Figure  7  -  The  shaded  region  shows  the  only  half  of  the  channel  that  needs  to  be 
analyzed  due  to  symmetry.  The  two  different  sets  of  axes  represent  the 
coordinate  systems  for  the  curved  section  and  the  straight  section,  whose 
solutions  need  to  be  spliced  together  at  the  intersection. 


For  the  curved  section  in  figure  7,  Laplace’s  equation  takes  the  following  form: 


v2o=o 


ld_ 

r  dr 


rJO 

a7 


+  - 


i  a2o 

r2  df 


=  0 


(7a) 

(7b) 


Assuming  that  <E> (r,(ft)  is  separable  into  independent  functions  of  r  and  (ft,  and  choosing  the 
solution  for  the  potential  that  has  decaying  exponentials  in  the  (ft  direction,  produces  a  potential  in  the 
curved  section  of  the  following  form: 


0(r,  <p)  =  [ +  A2e~k*  ] [A3  cos (k  In  r)  +  A4  sin (k  In  r )]  (8) 

Using  the  fact  that  the  electric  field  is  the  negative  gradient  of  the  potential  while  keeping  only 
decaying  exponentials,  and  using  the  boundary  conditions 

Er\r=a= 0  and  Er\r=b=0  (9) 

reduces  the  number  of  individual  coefficients  to  be  solved  and  results  in  the  following  series  solutions  for 
the  electric  field  components: 
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The  series  solutions  in  (10)  result  from  the  boundary  conditions  and  the  periodic  nature  of  the 
solution  basis  functions,  each  period  of  which  is  a  partial  solution  and  therefore  added  into  the  total  series 
solution.  Furthermore,  the  remaining  coefficients,  Cm  will  be  solved  when  the  curved  and  straight  section 
solutions  are  spliced  together. 

For  the  straight  section  seen  in  figure  7,  Laplace’s  equation  takes  the  following  form: 


v2o  =  o 

a2o  a2o 

-  +  - 


=  0 


(11a) 

(Hb) 


dx2  dy 2 

Again,  choosing  the  solutions  that  decay  from  the  intersection  line  leaves  us  with  the  potential  in 
the  straight  section: 

0(x,  y)  =  Bx  sin (kx)eky  +  B2  cos {kx)eky ,  (12) 

where  y  is  always  negative  as  defined  by  the  x-y  coordinate  system  in  figure  7,  so  the  exponents  in  (12) 
are  decaying  away  from  the  intersection  line. 

Taking  the  negative  gradient  of  (12)  and  using  that  to  create  a  Fourier  series  approximation  for 
the  correction  fields  provides  the  x  and  y  components  of  the  electric  field  as  follows: 


.  50  ^  ^  nK  .  ( nxx\  — 

y)  -  ~^r  =  — sin  —  r 

dx  „  w  V  w  ) 


rp  /  x  30  Y  T~v  n7l 

Ey  (x,  y )  =  =  2^~Dn  — cos 

dy  n  w 


/  \  nm 

'  nm ' 


V  W  J 


(13a) 


(13b) 


where  the  boundary  conditions 


Ex\x=n=0  and  Ex L=w=0 


(14) 
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result  in  the  series  solutions  and  reduce  the  number  of  coefficients,  leaving  only  Dm  which  can  be  solved, 
along  with  Cn  in  (10),  when  Ex,  Ey ,  Er ,  and  E $  are  spliced  at  the  intersection  line.  When  the  correction 
factor  components  in  (10)  and  (13)  are  superposed  with  the  abrupt-transition  1/r  field  solutions  from  (6) 
to  provide  the  correct  field  structure  far  from  the  intersection  line,  the  result  is  the  following  set  of 
equations: 
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(15) 


The  components  in  (15)  represent  the  field  in  both  the  straight  and  curved  sections,  and  are 
effectively  coupled  at  the  intersection  line  by  solving  for  the  remaining  coefficients,  Cn  and  Dn ,  using  the 
boundary  conditions  E^0=o=Ey\y=o  and  Er\^o=Ex\y=o-  The  number  of  coefficients  to  be  solved  will  depend 
on  the  number  of  terms  kept  in  each  of  the  series  in  (15).  Keeping  the  first  term  in  the  Fourier  series  for 
both  the  axial  and  transverse  components  in  the  straight  and  curved  sections  (represented  as  lxl,  where 
the  first  number  is  the  number  of  terms  kept  for  the  curved  channel’s  field  components  and  the  second 
number  is  the  number  of  terms  kept  for  the  straight  channel’s  field  components)  provides  sufficient 
accuracy  for  most  designs,  as  seen  in  Figures  8  and  9. 

To  find  the  coefficients  for  an  arbitrary  number  of  Fourier  components  in  each  section,  the  system 
of  equations  defining  the  coefficients  in  (15)  is  combined  in  a  relational  form  dependent  on  the 
appropriate  intersection  line  boundary  conditions  into  a  coupling  block  matrix.  This  creates  the  following 
linear  system,  used  to  solve  for  Cn  and  Dn\ 
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where  k,  are  simply  spatial  frequency  relations,  ifl/w,  m  is  the  number  of  terms  kept  in  the  series  for  the 
straight  section,  n  is  the  number  of  terms  kept  in  the  series  for  the  curved  section,  and  the  coefficients  £ 
77,  and  rare  solved  for  in  the  fairly  standard  method  of  determining  the  coefficients  of  Fourier  series. 


The  orthogonality  of  the  trigonometric  basis  functions  allows  a  particular  coefficient  to  be  isolated  by 
multiplying  the  entire  series  by  the  desired  basis  function  with  the  desired  frequency  and  integrating  over 
the  respective  period.  When  the  series  is  integrated  over  the  period,  it  provides  the  solution  for  the 
desired  coefficients,  £  V »  and  r,  which  are  related  to  each  other  by  the  linear  system  in  (16).  The 
generalized  results  for  the  coupling  matrix  coefficients  are: 
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Each  of  the  four  square  sub-blocks  in  the  coupling  matrix  of  (16)  has  dimensions  which  depend 
on  the  number  of  Fourier  terms  kept  in  each  of  the  curved  or  straight  sections.  The  two  left  sub-blocks 
are  mxn  matrices,  the  two  right  sub-blocks  are  mxm,  and  the  solution  vector  consists  of  ///+//  elements. 
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As  a  simple  example,  when  the  number  of  Fourier  components  is  equal  to  lxl,  the  following 


results  for  C\  and  D\  are  obtained: 
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Figure  8(a)  shows  the  magnitude  of  the  axial  electric  field  components  for  a  geometry  whose 
width-to-center  radius  ratio,  w/rc,  is  1.0,  which  is  at  one  extreme  end  of  the  feasible  design  space,  keeping 
one  Fourier  term  for  the  electric  fields  in  the  curved  and  straight  sections,  (lxl).  Figure  8(b)  shows  the 
results  of  keeping  three  Fourier  terms  for  the  electric  field  components  in  the  curved  and  straight  sections, 
(3x3).  Figure  9  shows  that  for  a  width-to-center  radius  ratio  of  0.4,  which  is  a  more  reasonable  upper 
limit  to  the  design  space  [4],  keeping  only  one  Fourier  component  achieves  excellent  agreement  with 
numeric  PDE  solutions  for  the  field  structure. 


1x1  FEM/Analytical  Plot  of  Component  Tangent 
to  Path  Length  from  Vo  to  O' 


(a) 


3x3  FEM/Analytical  Plot  of  Component  Tangent  to 
Path  Length  from  Vo  to  O' 


Figure  8  -  For  a  width-to-center  radius  ratio  ( w/rc )  of  1.0,  figure  8(a)  shows  the  axial  component  of  the  electric  field 
using  one  Fourier  term  for  curved  section  field  components  and  the  straight  section  field  components  (lxl).  figure 
8(b)  shows  the  same  comparison  using  three  Fourier  terms  for  the  curved  and  straight  section  field  components 
(3x3). 
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1x1  FEM/Analy  Plot  of  Component  'Tangent  to 
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3x3  FEM/Analy  Plot  of  Component  Tangent  to 
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Figure  9  -  For  a  width-to-center  radius  ratio  ( w/rc )  of  0.4,  figure  9(a)  shows  the  axial  field  component  result  for  lxl 
Fourier  terms,  and  figure  9(b)  shows  the  results  for  3x3  Fourier  terms.  For  smaller  w/rc  ratios  the  accuracy  of  the 
analytic  field  expressions  in  (15)  increases  for  a  given  number  of  Fourier  terms. 


One  final  observation  that  can  be  made  from  the  Laplacian  derivation  is  the  persistence  of  the 
transition  field  components  away  from  the  intersection  line  of  the  channel  sections.  This  persistence  is 
important  for  measuring  the  turn-induced  variance  because  the  band  must  be  outside  of  the  transition 
fields  before  the  turn-induced  variance  is  complete  and  a  measurement  can  be  made.  Looking  at  (15),  it 
can  be  seen  that  the  transition  components  have  exponential  decay  factors,  which  can  be  used  to 
determine  how  far  the  transition  fields  persist  from  the  intersection  line.  The  transition  fields  fall  to  4% 
of  their  original  value  at  about  a  channel’s  width  length  into  the  straight  section,  and  they  fall  to  about  4% 
of  their  original  value  at  about  a  channel’s  width  arc-length  into  the  curved  section  along  the  center 
radius.  This  decay  rate  is  for  the  first  harmonic  of  the  Fourier  series,  which  is  the  least  attenuated  of  all 
the  harmonics.  So  after  a  channel’s  width  length  in  both  directions  from  the  intersection  line  along  the 
center  radius  (shown  in  figure  10),  the  transition  fields  have  decayed  so  that  the  field  structure  is 
essentially  defined  by  the  abrupt-transition  1/r  field  description. 
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Figure  10  -  The  persistence  of  the  transition  fields  from  the  intersection 
line,  highlighted  in  the  green  boxes,  for  geometry  with  a  width-to-center 
radius,  w/rc  ratio  of  1.0  on  the  left  and  a  w/rc  of  0.4  on  the  right. 

To  make  the  final  comparison  with  the  abrupt-transition  field  results  seen  in  figure  6,  figure  11 
shows  the  axial  and  transverse  components  of  the  Laplacian  field  structure,  for  w/rc=  1.0  along  axial  lines 
at  w/ 4,  wl 2,  and  3w/4,  which  now  accounts  for  the  missing  field  structure. 

In  summary,  this  chapter  has  shown  two  models  of  the  electric  field  structure  for  serpentine 
channels  gradually  reducing  the  assumptions  found  in  [2, 3, 4, 5].  The  effect  the  new  field  structure  has  on 
the  variance  of  a  band  will  be  shown  in  the  next  chapter. 


Laplacian  (1x1):  Analytical  Plot  of  Component  Tangent  to 
Path  Length  from  Vo  to  O' 


Laplacian  (1x1):  Analytical  Plot  of  Component  'Normal  to 
Path  Length  from  Vo  to  O' 


Figure  1 1  -  A  plot  of  the  axial  (on  the  left)  and  transverse  (on  the  right)  components  of  the  electric  field  as  described  by 
the  Laplacian  derivation  in  Eq.  8.  Compared  to  Figure  2,  there  is  a  smoother  transition  between  the  straight  and  curved 
section  field  structure  on  the  left,  and  the  transverse  components  on  the  right,  missing  in  Figure  2,  are  now  present. 
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3.  Analyte  Dispersion  Results 


This  chapter  will  use  the  electric  field  models  developed  in  the  previous  chapter  to  create  models 
of  the  variance  of  the  analyte  band  as  it  travels  through  the  serpentine  separation  channel,  which  is  a 
common  separation  topology  [3].  The  serpentine  channels  are  simply  complementary  connections  of  the 
single  turns  and  straight  sections  seen  in  the  previous  chapter,  as  shown  in  figure  12.  The  method  of 
measuring  the  variance  of  the  band  will  be  laid  out,  followed  by  a  description  of  the  various  operating 
regimes  of  the  separation  system.  A  variance  model  using  the  abrupt-transition  1/r  fields  developed  in  the 
previous  chapter  will  then  be  created  and  used  to  generate  a  model  using  the  continuous-transition 
Laplacian  fields,  also  developed  in  the  previous  chapter.  Comparisons  to  FEMLAB  numerical  PDE 
simulations  are  made  as  they  are  needed.  Finally,  there  are  several  sources  of  dispersion  in  serpentine 
separation  systems,  such  as  the  injector,  detector,  electric  field  structure,  geometry  and  diffusion.  The 
models  in  this  chapter  focus  on  the  electric  field  structure  and  geometry  as  sources  of  dispersion,  since 
such  sources  arise  from  serpentine  topologies. 

3.1  Measuring  Dispersion 

The  process  of  measuring  the  dispersion  of  a  continuous  concentration  distribution  is  essential  to 
creating  models  of  separation  systems.  This  section  introduces  the  measurement  method  used  for  the 
models  in  this  report. 

Real  systems  measure  the  dispersion  of  an  analyte  band  in  a  channel  with  various  sensors  ranging 
from  2-dimensional  (2-D)  CCD  arrays  to  laser  scanning  single-point  optical  devices  [12].  These  sensor 
devices  have  practical  experimental  limitations  such  as  the  resolution  of  the  pixels  in  the  CCD  and  the 
wavelength  of  the  laser-scanning  element.  These  limitations  are  usually  characterized  and  modeled  as  an 
additional  source  of  dispersion  in  the  system.  Since  injectors  and  detectors  are  not  the  focus  of  this  work, 
only  the  effects  of  the  serpentine  channels  are  included  in  the  calculation  and  derivation  of  models  of  the 
band  dispersion.  Injection  and  detector  variance  models  will  be  the  focus  of  later  work. 
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For  a  continuous  2-D  concentration  distribution,  the  method  of  determining  the  analyte  variance 
is  summarized  in  figure  12.  The  concentration  distribution  is  reduced  to  a  y-dependent  1-dimensional  (1- 
D)  concentration  profile  by  averaging  along  the  x-axis.  Then  the  variance  of  this  1-D  profile  is  found 
using  the  standard  formula  shown  in  figure  12. 


<j2  = 


 J(y-y0): 


cavg(y)dy 


cavg  (. y)dy 


Figure  12  -  The  process  of  converting  a  concentration  distribution  into  a  measure  of  variance. 


3.2  Non-Dimensional  Serpentine  Channel  Parameters 

This  section  will  first  show  a  non-dimensionalization  of  the  advection-diffusion  equation,  which 
describes  the  dynamics  of  the  separation  channels.  Then  a  complete  set  of  non-dimensional  characteristic 
parameters  will  be  derived  using  the  Buckingham  71-Theorem.  These  steps  are  useful  for  gaining  an 
intuitive  understanding  of  how  the  separation  systems  operate. 

The  dynamics  of  the  analyte  band  traveling  through  a  channel  are  described  by  the  advection- 
diffusion  partial  differential  equation,  repeated  from  (4)  here: 

—  +  u  •  Vc  =  DW2c  (19a) 

dt 

Non-dimensionalizing  (19a)  gives 

dc  7  ^  *  1  * — 7  =^2^; 

—  +  u  -V  c  =— V  c  ,  (19b) 

dt  Pe 

where 
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In  the  non-dimensional  form  of  (19b),  the  balance  between  the  advection  and  diffusion  terms 
more  evidently  describes  the  dynamics  of  a  system.  In  this  case,  the  Peclet  number  describes  this  balance 
between  advection  and  diffusion,  where 


Pe  = 


u0w 

D  ’ 


(21) 


and  the  concentration,  length,  and  time  dimensions  have  been  scaled  by  characteristic  parameters,  shown 
in  (20),  such  as  the  initial  concentration,  cQ,  the  width  of  the  channel,  w,  and  the  ratio  of  the  channel  width 
to  mean  velocity,  w/u0,  respectively.  Here  the  mean  velocity,  uoy  is  the  mobility  times  the  abrupt- 
transition  1/r  field  value  in  the  straight  sections,  Ey,  which  is  the  mean  field  strength  for  the  entire 
channel,  including  the  curves.  This  type  of  dimensional  analysis  in  (19)  can  be  used  to  identify  four  basic 
regimes  in  which  the  system  can  operate,  depending  on  its  characteristic  non-dimensionalized  parameters. 
In  this  case,  the  Peclet  number  is  one  of  the  characteristic  non-dimensional  parameters. 

To  discover  any  additional  non-dimensional  parameters  that  characterize  the  operating  regime, 
the  Buckingham  71-Theorem  can  be  utilized,  as  shown  in  the  Appendix.  Based  on  the  physical  parameters 
describing  the  band  variance  of  a  serpentine  microfluidic  system,  such  as  the  channel  width  (w),  turn 
center  radius  (rc),  length  between  complementary  turns  (Lbt),  the  electric  field  strength  ( E0 ),  diffusivity 
( D ),  and  species  mobility  (//),  the  following  maximal  set  of  three  non-dimensional  terms  fully 
characterize  the  system:  rc/w,  LJw ,  and  u0w/D.  Where  u0w/D  is  recognized  as  the  Peclet  number  based 
on  channel  width.  Defining  the  non-dimensional  parameters  Tbt=wPe/Lbt  (the  ratio  of  the  time  it  takes  a 
particle  to  travel  the  width  of  the  channel  to  the  time  it  takes  a  particle  to  travel  the  length  between  turns) 
and  Tt=wPel(7rrc)  (the  ratio  of  the  time  it  takes  a  particle  to  travel  the  width  of  the  channel  to  the  time  it 
takes  a  particle  to  travel  the  length  of  the  turn  along  rc),  three  equivalent  non-dimensional  parameters  are: 
Tbh  Th  and  w/rc.  These  parameters  can  be  used  to  fully  characterize  the  regimes  in  the  problem  space  of 
the  separation  system,  as  seen  in  the  next  section. 
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3.3  System  Dispersion  Regimes 

This  section  will  introduce  the  various  regimes  in  which  the  separation  system  can  operate, 
followed  by  several  examples  using  FEMLAB. 


Figure  13  -  FEMLAB  simulations  created  that  span  the  Thl  space  crossing  four 
identifiable  regimes,  which  go  from  left  to  right,  red  to  blue,  as  the  pure 
diffusion,  advection-diffusion  transition,  pure  advection,  and  high-pure 
advection  regimes. 

The  four  regimes  are  easily  identified  in  the  graph  by  the  three  knees  separating  the  regimes  with 
four  different  characteristic  slopes.  The  four  regimes  from  left  to  right  are:  Pure  diffusion,  advection- 
diffusion  transition,  pure  advection,  and  high-pure  advection.  The  pure  diffusion  regime  is  characterized 
by  the  band  being  always  homogenous  orthogonal  to  the  channel  axis  and  widening  along  the  channel 
axis.  The  advection-diffusion  transition  region  is  characterized  by  the  band  being  mostly  homogenous, 
but  experiencing  some  skewing  from  advection  around  the  turns.  However,  in  this  regime  there  is  not 
enough  advection  dominance  to  permit  unskewing  after  complementary  turns.  In  the  pure  advection 
regime,  the  band  strongly  expresses  skewing  and  advection  is  now  dominant  enough  to  permit  some 
unskewing  of  the  band,  i.e.,  some  dispersion  reduction,  however,  axial  band  widening  is  still  present.  In 
the  high-pure  advection  regime,  the  analyte  bands  strongly  express  skewing  and  unskewing  after  turns,  a 
characteristic  parabolic-like  shape  after  complementary  turns,  and  no  significant  widening  of  the  band 
from  diffusion.  Figures  14a  through  14d  show  examples  of  bands  traveling  in  each  of  these  regimes. 
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Figure  14a  -  Examples  of  the  various  regimes.  This  figure  shows  the  pure  diffusion  regime  in  the  red  section  of 
figure  13.  The  band  shows  no  skewing  after  turns  and  maintains  essentially  a  1-D  concentration  gradient  along  the 
channel  axis.  Note  the  band  has  traveled  from  the  middle  of  the  first  leg  to  the  middle  of  the  last  leg. 


Table  1  -  System  parameters  for  figure  14a. 
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Figure  14b  -  Examples  of  the  various  regimes.  This  figure  shows  the  advection-diffusion  transition  regime  in  the 
yellow  section  of  figure  13.  There  is  skewing  after  a  turn,  but  significant  transverse  diffusion  remains  that  keeps  the 
band  fairly  homogeneous  across  the  channel  width  and  significant  axial  diffusion  spreads  the  band  axially. 
Therefore,  complementary  turns  cannot  unskew  the  band,  but  instead  “overskew”  the  band.  Note  the  band  has 
traveled  from  the  middle  of  the  first  leg  to  the  middle  of  the  last  leg. 


Table  2  -  System  parameters  for  figure  14b. 
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Figure  14c  -  Examples  of  the  various  regimes.  In  the  pure-advection  regime  shown  here,  representing  the  green 
region  of  figure  13,  there  is  significant  skewing  of  the  band  after  a  turn  and  transverse  diffusion  is  small  enough  to 
prevent  the  band  from  homogenizing  across  the  channel  width.  The  field  effects  are  dominant  enough  over 
diffusion  to  allow  some  unskewing  after  complementary  turns,  some  “over skewing”  is  still  present,  and  axial 
diffusion  is  still  present  enough  to  spread  the  band  axially.  Note  the  band  has  traveled  from  the  middle  of  the  first 
leg  to  the  middle  of  the  last  leg. 


Table  3  -  System  parameters  for  figure  14c. 
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Figure  14d  -  Examples  of  the  various  regimes.  The  last  regime,  shown  in  this  figure  and  representing  the  blue 
region  in  figure  13,  is  the  high-pure  advection  regime.  In  this  regime,  the  effects  of  diffusion  are  negligible,  so  there 
is  no  significant  widening  of  the  band  axially  or  transversely.  The  effects  of  the  geometry  and  field  structure  define 
the  skew  and  unskew  of  the  band.  The  unskewing  of  the  band  after  complementary  turns  results  in  the  characteristic 
parabolic  shape.  Note  the  band  has  traveled  from  the  middle  of  the  first  leg  to  the  middle  of  the  last  leg. 


Table  4  -  System  parameters  for  figure  14d. 
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Most  microfluidic  designs  lie  within  the  pure  and  high-pure  advection  regimes  (green  and  blue) 
[3].  Out  of  these  two  regimes,  the  majority  of  systems  would  fall  into  the  pure-advection  regime. 
However,  there  are  still  some  very  important  examples  that  fall  into  the  high-pure  advection  regime.  For 
example,  a  typical  example  of  a  system  found  in  the  high-pure  advection  regime  is  DNA  electrophoresis 
where  a  sieving  matrix  is  employed,  as  seen  in  [12].  A  typical  Peclet  number  for  such  systems  is  around 
15,000.  Within  the  scope  of  these  regimes,  the  next  section  will  present  analysis  starting  with  the  abrupt- 
transition  Hr  field  structure  to  analytically  determine  the  variance  of  a  band  as  it  travels  through 
serpentine  channels  within  the  high-pure  advection  regime.  Further  sections  will  begin  to  look  at  the 
regimes  with  where  diffusion  is  more  important. 

3.4  High-Pure  Advection  Analytical  Dispersion  with  1 -Assumption  Field  Model 

Using  the  concepts  introduced  in  the  previous  sections  of  this  chapter,  this  section  will  produce  an 
analytical  model  of  the  high  Peclet  dispersion  in  serpentine  channels  with  any  number  of  turns  and  the 
abrupt-transition  Hr  electric  field  structure.  This  model  will  then  be  used  to  explore  several  previously 
unpredicted  physical  phenomena. 

For  bands  traveling  though  the  serpentine  channels  in  the  high-pure  advection  regime,  we  can 
ignore  the  effects  of  diffusion.  Band  variance  is  therefore  determined  only  by  the  channel  geometry  and 
field  structure.  To  analyze  the  abrupt-transition  field  dispersion,  only  the  curved  sections  need  to  be 
considered,  since  the  uniform  fields  in  the  straight  sections  do  not  contribute  to  the  dispersion  when 
diffusion  is  ignored.  To  find  the  variance  of  a  band  traveling  through  a  serpentine  channel,  the  band  can 
be  thought  of  as  individual  ‘particles’  following  an  axial  path  through  the  channel.  The  differences  in  the 
travel  time  through  the  curves  for  ‘particles’  spread  across  the  width  of  the  channel  will  serve  as  a 
measure  of  the  band  variance  after  traveling  through  the  channel. 
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Figure  15  -  The  difference  in  times  of  outer  and  inner  particles  traveling  through  a  single  turn. 


The  time  it  takes  a  particle  at  radius  r  to  travel  through  the  channel  is  the  distance  divided  by  the 
velocity, 


T(r)  = 


7ir 

ME4,(r)  ’ 


(22) 


where  fi  is  the  particle  mobility  and  E $  is  the  circumferential  field  in  the  curve  using  the  abrupt-transition 
description  in  (6b),  which  also  depends  on  r.  So  the  particle  along  the  inner  most  radius,  a ,  takes 
Ti=m/[juE^a)]  seconds  and  the  particle  along  the  outer  most  radius  takes  TQ=7ibl[juE(/^b)]  seconds.  Their 
difference  is  To-T^TiZ/Aib/E^bya/E^a)],  which  is  the  time  spent  by  the  band  traveling  in  the  straight 
section  immediately  after  the  curve.  Therefore,  the  time  any  given  particle  spends  in  the  straight  section 
after  the  curve  is  given  by: 


A  T(r)  = 


71 


b 


EAb )  EAr) 


(23) 


This  AT,  in  (23),  multiplied  by  the  uniform  velocity  in  the  straight  section  gives  the  distance  traveled  by 
each  particle  in  the  uniform  straight  section: 


AUn  =  L,EATU)  = 


xln(bla)[b2_r2h 


w 


(24) 


where  a ,  b ,  and  w  are  the  dimensions  as  defined  in  figure  5b. 

Using  (24),  the  variance  of  the  profile  about  its  mean  using  the  sample  variance, 
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(25) 


A  similar  procedure  as  above  can  be  used  to  derive  the  variance  for  any  number  of  turns.  For  this 
derivation,  however,  the  absolute  time  spent  in  the  curve,  (22),  and  its  effect  on  the  displacement  in  the 
straight  section  after  the  turns  is  used  to  calculate  the  variance.  Since  the  variance  is  a  centered  moment, 
whether  or  not  the  outer  radius  length  or  time  is  subtracted  out  has  no  effect  on  the  resulting  variance. 
This  greatly  simplifies  the  algebra,  but  does  not  give  physical  significance  to  the  real  displacement  of  the 
particles,  since  using  the  time  spent  in  the  curves  will  result  in  a  displacement  in  the  straight  section 
appearing  inverted.  For  example,  for  a  single  turn,  an  inner  radius  particle  has  the  shortest  travel  time 
through  the  turn,  so  it  will  have  the  shortest  length  afterwards  if  this  time  spent  in  the  curves  is  used  to 
extend  the  particles  into  the  straight  section.  This  is  obviously  inverted  from  the  physical  result,  since  the 
innermost  particle  should  be  the  furthest  ahead  in  the  straight  section.  So,  a  physical  interpretation  of  the 
particle  displacements  requires  that  the  time  used  to  extend  the  particles  into  the  straight  section  is  the 
outer  radius  particle  time  minus  the  particle  time  traveling  at  a  desired  radius,  as  mentioned  for  the  results 
in  (24)-(26).  For  the  general  case  of  n  turns,  the  physical  displacements  of  the  particles  is  given  by: 

A L{r)  =  - ln^ - ^  K  (b2  -  r 2 )  +  N2  ((w  +  2a  - b)2  -(w  +  2a-r)2)]  (27) 

w 

However,  to  find  the  variance  using  a  simpler  form  of  displacement,  the  time  spent  in  the  curves 
is  utilized  to  find  the  distance  traveled  in  the  straight  section,  which  results  in  the  following: 

L(r)  =  ^ln(^/a)  [jV)r2  +  N2(w+2a  —  r)2\,  (28) 

w 

Since  (28)  is  simply  a  constant  minus  (27),  they  both  have  the  same  variance.  (28)  is  then  used  in  the 
sample  variance  equation: 
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(29) 


a 


which  gives  the  following  result  for  the  variance  with  the  abrupt-transition  1/r  electric  fields  for  n 
complementary  turns: 


2  n  In  (b/a) 

& ab  ~  3 
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D 

J  N,r2  +N2(w  +  2a-r)2  - 


3  \  A 


(N,  +N2)((w  +  ay  -a") 
3  w 


dr , 


(30) 


where 


v  = 


n  +  mod(n,2) 


N2  = 


n  +  mod  (n  + 1,2)  - 1 


(31) 

(32) 


and  r  is  the  radius  starting  in  the  first  curve,  n  is  the  number  of  turns,  and  mod(n,2)  is  the  modulus  after 
division  function  that  gives  the  remainder  of  nil.  To  see  the  complete  expanded  solution  of  (30),  refer  to 
the  Appendix. 

A  quick  example  of  (30)  shows  that  for  a  single  set  of  complementary  turns  the  following  fairly 
simple  result  is  produced: 

(y]b  =^7t2w2\n2(bla)  (33) 

Several  interesting  results  are  now  apparent.  First,  figure  16  shows  the  resulting  band  edge  after  a 
single  set  of  complementary  turns.  It  shows  that  with  the  1 -assumption  electric  field  model  with  the 
abrupt  transition  in  the  field  structure  between  channel  sections,  unskewing  is  not  complete  and  leaves  a 
parabolic-like  shape  after  complementary  turns.  Whereas  the  2-assumption  models  in  [2, 3, 4, 5]  with  the 
abrupt  transition  in  the  field  between  channel  sections  and  the  linear  electric  field  in  the  turn,  predict  a 
complete  unskewing  of  the  band  after  complementary  turns. 
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Figure  16  -  Band  shape  after  a  set  of  complementary  turns,  traveling 
from  left  to  right.  A  distorted  band  results,  unlike  when  using  a  linear 
description  of  the  field  structure.  A  Matlab  particle  simulation  created 
these  results,  which  will  be  fully  described  in  the  next  section. 

The  resultant  band  shape  is  symmetric  about  the  center  of  the  channel,  but  part  of  what 
contributes  to  the  parabolic  shape  is  the  non-symmetry  of  the  electric  field  in  the  turn  about  the  central 
channel  axis  and  such  interaction  from  turn  to  turn.  This  non-symmetry  occurs  because  the  field’s  mean 
value  is  not  precisely  in  the  center  of  the  channel,  but  at  a  radius  less  than  the  center.  Since  the 
complementary  turns  mirror  the  curvature  of  the  previous  turn,  but  the  fields  are  different  at  the  center  of 
the  channel,  the  center  of  the  band  gets  and  stays  further  ahead  than  the  two  band  edges,  which  both 
symmetrically  experience  the  long  path-length  and  lessened  field-strength  of  the  outer  radius  of  the 
channel.  If  the  field  were  symmetric  about  the  center  of  the  channel,  then  each  particle  would  have 
experienced  the  same  electric  field  as  its  mirrored  partner  from  the  previous  channel,  thereby  equalizing 
the  effects  from  the  previous  turn  and  flattening  the  band.  This  is  the  case  with  any  linear  velocity 
profile,  as  opposed  to  any  non-linear  profile  that  is  not  symmetric  about  its  mean. 

Another  interesting  observation,  shown  in  figure  17,  is  that  as  the  number  of  complementary  turns 
continues  to  grow,  the  resulting  variance  from  the  parabolic-like  band  shape  after  each  set  of 
complementary  turns  grows  quadratically.  This  is  for  the  high-pure  advection  regime,  where  the  effects 
of  diffusion  are  negligible. 
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Figure  17  -  Variance  growth  with  increasing  number 
of  complementary  turns  in  system.  w/rc=0.4 ,  a=10pm. 


Figure  18  -  Plot  of  the  turn  dispersion  after  even  and  odd 
turns  for  w/rc=0A.  The  top  line  represents  the  dispersion  after 
only  the  odd  numbered  turns.  The  bottom  line  represents  the 
dispersion  after  only  the  even  numbered  turns.  The  zig¬ 
zagging  line  represents  the  variance  after  both  sets  of  turns. 


The  source  of  this  quadratic  growth  is  apparent  in  the  expanded  form  of  (30),  seen  in  (al9)  in  the 
Appendix.  Every  term  of  (al9)  contains  an  n2,  which  indicates  that  the  variance  grows  quadratically  with 
both  even-  and  odd-numbered  turns.  An  example  of  this  is  shown  in  figure  18. 

The  variance  grows  quadratically  for  either  even  or  odd  sets  of  turns.  However,  when  following 
the  variance  after  both  sets  of  turns,  it  alternates  high  and  low,  since  the  variance  grows  after  odd  turns 
and  is  subsequently  reduced  after  even-complementary  turns.  Figure  18  makes  it  apparent  why  odd- 
turned  designs  are  not  practical,  since  they  always  produce  more  variance  than  an  even-turned  design. 

In  reality,  since  the  system  is  limited  in  size,  the  variance  cannot  continue  to  grow  to  an 
arbitrarily  large  quantity.  However,  until  a  physical  limit  is  reached,  the  growth  rate,  since  quadratic,  will 
be  significantly  large.  Such  effects  will  be  discussed  further  in  section  3.6.  The  variance  results  in 
figures  17  and  18  are  all  valid  for  high  Peclet  regimes  (the  high-pure  advection  region  in  figure  13).  In 
the  next  section,  diffusion  will  be  considered  to  examine  the  pure  advection  (green)  regime  of  figure  13, 
where  many  microfluidic  designs  lie  [3,4].  Even  with  diffusion,  super-linear  growth  rates  for  dispersion 
can  be  observed. 
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3.5  Particle  Simulation  Results  with  0-Assumption  Field  Model 

With  the  full  Laplacian  field  structure,  a  direct  first-principles  derivation  of  a  closed-form 
algebraic  expression  for  the  high-Peclet  band  dispersion  is  not  possible  using  the  analysis  from  the 
previous  section  due  to  the  non-fully  developed  flow  patterns  created  by  the  transition  electric  fields.  So 
in  this  section,  a  particle  simulator  was  created  in  Matlab,  [13],  using  numerical  solutions  of  the  full 
Laplacian  field  equations,  (15),  to  determine  the  dynamics  of  particles  placed  in  the  channel.  These 
dynamics  are  described  by  the  Helmholtz-Smoluchowski  equation,  re-written  from  (3)  here: 

dx 

—  =  //£,  (34) 

dt 

where  vis  the  position  vector,  //  is  the  species  mobility,  and  E  is  the  electric  field,  which  can  be 
represented  by  the  full  Laplacian  solution,  (15),  or  the  abrupt-transition  solution,  (6). 

The  particle  simulator  is  fed  the  geometry  and  physical  properties  of  the  system  such  as  the  inner 
and  outer  radii,  a  and  b ,  the  straight  leg  lengths,  L,  the  number  of  curves,  n ,  the  particle  mobility,  //,  the 
applied  electric  potential,  V0,  and  the  number  of  particles,  pts.  Then,  the  simulator  numerically  solves 
(34)  to  provide  the  particle  trajectories.  Figure  19  shows  the  geometry  input  for  the  simulator. 


Figure  19  -  Geometry  input  to  the  particle  simulator. 

After  running  the  time  necessary  for  the  mean  of  the  band  to  reach  the  defined  end  of  the 
serpentine  channel,  the  particle  simulator  produces  an  array  of  particle  locations.  Then,  by  using  the 
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discretized  form  of  (25),  implemented  through  the  Matlab  var( )  function,  the  sample  variance  is 
calculated,  which  represents  the  variance  of  the  band  in  the  channel. 

Now,  the  particle  simulator  will  be  used  to  verify  the  analytical  abrupt-transition  1/r  field 
variance.  Then,  simulation  results  can  be  generated  for  the  full  Laplacian  transition  field,  and  new 
comparisons  can  be  made  between  the  abrupt  and  Laplacian  transition  field  variance  predictions. 


Figure  20  -  (a)  Comparison  of  the  abrupt- transition  field  particle  simulation  variance  to  the  analytical  expression  of 
the  variance  for  a  single  turn  with  w/rc  ranging  from  0.1  to  1.0.  (b)  Comparison  of  the  particle  abrupt-transition  field 
particle  simulation  variance  to  the  analytical  expression  of  the  variance  for  complementary  turns  with  w/rc  from  0.1 
to  1.0.  In  both  cases,  agreement  is  excellent. 

Figures  20a  and  20b  compare  and  verify  the  particle  simulator’s  abrupt-transition  Hr  field 
variance  results  to  the  analytical  expressions  from  the  previous  section.  Since  w/rc  is  the  only  non- 
dimensional  variable  for  infinite-Peclet  dispersion,  which  can  be  verified  with  the  Buckingham  71- 
Theorem  in  the  Appendix,  the  first  graph  shows  the  results  for  a  range  of  w/rc  values  from  0.1  to  1.0, 
which  more  than  encompasses  the  majority  of  the  practical  design  space  which  lies  in  the  range  0.1  to  0.4, 
as  previously  mentioned.  From  these  figures,  the  analytical  expressions  from  the  previous  section  can  be 
seen  to  be  in  excellent  agreement  with  the  particle  simulator  across  the  entire  range  of  w/rc,  with  the  error 
consistently  near  the  floor  of  the  computer’s  precision. 

Next,  dispersion  results  for  the  full  Laplacian  transition  fields  can  be  generated  for  the  same  range 
of  w/rc  values  (0.1  to  1.0)  and  compared  to  the  abrupt-transition  field  dispersion.  Figure  21  shows  the 
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results  of  a  particle  simulation  at  a  single  w/rc  value,  0.667,  which  runs  with  abrupt-transition  (red)  and 
full  Laplacian  transition  (blue)  electric  fields.  They  are  compared  with  the  variance  results  from  a  linear 
(green)  electric  field  in  the  curve,  which  comes  from  the  expression  skew2/ 12. 


Figure  21  -  Particle  simulations  traveling  left  to  right  for  w/rc=0.667  with 
the  full  Laplacian  transition  field  (blue)  and  the  abrupt-transition  1/r  field 
(red)  compared  with  the  linear  field  approximation  skew  (green). 

Figure  22a  compares  the  abrupt-transition  1/r  and  Laplacian  field  variance  results  with  the  skew 
variance  resulting  from  the  2-assumption  linear  field  model  for  w/rc  ranging  from  0.1  to  1.0  for  a  single 
turn.  Figure  22b  shows  the  two  results  compared  with  the  2-assumption  linear  field  model  for  a  set  of 
complementary  turns  for  the  same  wtrc  range.  In  figure  22a,  the  difference  between  the  two  descriptions 
is  not  that  significant  or  that  different  from  the  linear  field  variance,  skew2/l2 ,  for  small  w/rc.  This  region 
of  w/rc  is  where  the  error  of  the  linear  approximation  for  the  field  in  the  curve  is  becoming  smaller,  so  the 
agreement  of  the  three  becomes  fairly  good.  However,  for  complementary  turns,  even  at  low  w/rc  the 
linear  and  1/r  abrupt-transition  approximations  miss  vital  contributions  to  the  dispersion  and  therefore  do 
not  agree  well  with  the  full  Laplacian  transition  electric  field  model. 
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Figure  22  -  The  figure  on  the  left  shows  the  resulting  variance  for  single  turn  channels  with  a  w/rc  from  0. 1  to 
1.0.  The  figure  on  the  right  shows  the  resulting  variance  for  complementary  turns  for  the  same  range  of  w/rc. 

In  the  previous  section  with  abrupt-transition  1/r  fields,  a  quadratic  growth  was  discovered  in  the 
variance  after  traveling  through  n  complementary  sets  of  turns  using  the  analytical  equations  for  abrupt- 
transition  field  dispersion.  Using  the  particle  simulator,  the  same  test  can  be  run  using  the  Laplacian 
transition  fields.  Figure  23a  shows  the  results  of  a  particle  simulation  for  the  abrupt-transition  1/r  and 
Laplacian  fields  through  three  sets  of  complementary  turns.  The  increase  in  band  variance  is  apparent 
after  each  individual  turn  and  particularly  after  each  set  of  complementary  turns.  Even  after  the  odd  turns 
that  come  after  at  least  one  set  of  complementary  turns,  the  mid-section  of  the  band  begins  to  outrun  the 
edges,  which  occurs  when  the  mid-section  gets  far  enough  ahead  of  the  band  edges  so  that  not  even  the 
larger  velocity  of  the  inner  radius  particles  is  enough  to  catch  up  with  the  forward  stretched  mid-band 
particles.  Figure  23b  shows  the  results  of  this  test  for  six  complementary  sets  of  turns,  or  twelve  turns 
total.  Again,  quadratic  growth  is  apparent.  Of  course,  the  linear  field  approximations  would  predict 
exactly  zero  variance  after  each  set  of  complementary  turns. 
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Complementary  Turn  Dispersion 


(b) 


Figure  23  -  (a)  Three  complementary  turns  simulated  with  the  particle  simulator  for  w/rc=0A.  (b)  Variance 
after  multiple  complementary  turns  for  Laplacian  and  abrupt- transition  1/r  fields,  showing  quadratic  growth 
in  both  cases.  w/rc=0A,  a=  lOqm. 


One  observation  from  figure  23b  is  that  since  both  Laplacian  and  abrupt-transition  Hr  field 
variance  descriptions  have  quadratic  growth  in  this  high-pure  advection  regime,  they  differ  only  by  a 
constant,  which  in  the  case  of  figure  23b  is  slightly  less  than  four.  This  fact  will  be  used  in  section  3.7  to 
determine  a  relationship  between  analytical  variance  descriptions  using  the  abrupt-transition  electric  field 
in  (30)  and  the  variance  descriptions  using  the  Laplacian  electric  field  determined  by  the  particle 
simulator  in  order  to  create  a  closed  form  analytical  description  of  the  Laplacian  variance  in  the  high-pure 
advection  regime. 

To  validate  some  of  the  Laplacian  particle  simulation  results  and  to  extend  them  to  include 
diffusion  in  the  other  three  regimes,  FEMLAB  simulations  will  be  shown  in  the  next  section. 


3.6  Multi-Regime  FEMLAB  Validation  Results 

FEMLAB  can  be  used  to  numerically  solve  the  partial  differential  advection-diffusion  equation 
and  the  electrostatic  field  structure  necessary  to  solve  the  dynamics  of  a  separation  system  [14].  In  this 
section,  FEMLAB  will  be  used  to  provide  validation  results  for  the  models  shown  in  previous  sections 
and  self- validation  compared  to  published  experiments.  The  FEMLAB  simulations  involve  an 
unsymmetric  coupling,  where  the  advection  part  of  the  advection-diffusion  equation,  (4),  is  coupled  to  the 
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electric  field  through  the  Helmholtz-Schmolochowski  relation,  (3),  but  according  to  the  electroneutrality 
condition  previously  assumed,  the  species  concentration  has  no  effect  on  the  electric  field.  FEMLAB  can 
also  be  used  to  couple  any  number  of  custom-input  equations.  Using  this  feature,  FEMLAB  can  be  used 
to  verify  the  abrupt-transition  field  results  as  well. 

The  first  simulation  run  with  FEMLAB  is  a  complementary  turn,  with  w/rc=0.4,  in  the  high-pure 
advection  regime  that  compares  to  the  particle  simulations  for  abrupt-transition  and  Laplacian  fields. 

Abrupt-Transition  Hr 


Figure  24  -  FEMLAB  validation  of  abrupt-transition  and  Laplacian 
field  particle  simulations,  where  FEMLAB  used  a  diffusion 
coefficient  of  zero.  This  simulation  uses  a  center  radius,  rc-  12.5pm, 
length  between  turns,  L^=80pm,  mobility,  //=1.5-10‘7m2/Vs,  straight 
leg  electric  field,  £^=8420.1 7 V/m,  and  diffusivity,  D- 0. 

Figure  24  shows  this  validation  by  using  a  diffusion  coefficient  of  zero  to  obtain  an  infinite  Peclet 
number,  or  as  close  to  infinite  as  the  solver  tolerances  allow.  The  FEMLAB  simulation  confirms  the 
different  shapes  after  a  set  of  complementary  turns,  indicating  that  the  Laplacian  fields  provide  a 
significant  increase  in  variance  over  the  simpler  abrupt-transition  fields,  and  an  infinite  increase  compared 
to  the  linear  fields,  which  would  completely  unskew  the  band.  The  spreading  of  the  band  apparent  in  the 
FEMLAB  figures  is  a  result  of  the  initial  Gaussian  band  that  was  injected  into  the  channels,  in  contrast  to 
the  uniform  initial  band  that  is  implicit  in  the  particle  simulation. 
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Now  to  compare  FEMLAB  with  experimental  results,  a  simulation  can  be  run  that  replicates  the 
system  in  Ramsey’s  experiment  in  [4].  This  involves  measuring  the  dispersion  as  the  band  moves 
through  the  channel,  excluding  the  curved  sections.  The  system  has  an  rc=125|am,  a  w=60pm,  an 
Lbt=  6mm,  an  average  electric  field  intensity  of  60V/cm,  a  mobility  of  1 .195*  10'8  m2/Vs,  and  a  Peclet 
number  of  138.  The  geometry  is  shown  in  figure  25a,  and  the  results  of  the  simulation  are  shown  in 
figure  25b. 


I _ Lbt= 6mm _ ^ 

*  <^== 

w-60gm 


(a) 


Variance  Versus  Time 
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Figure  25  -  (a)  The  geometry  corresponding  to  Ramsey’s  experiment  in  [4].  (b) 
FEMLAB  simulation  of  a  separation  based  on  Ramsey’s  experiment  using  the  geometry 
in  figure  25a  for  a  Pe=138. 


The  slope  of  the  lines  in  figure  25b  corresponds  to  2D,  where  D  is  the  diffusion  coefficient  and 
2D  is  the  rate  of  variance  growth,  and  the  jumps  after  the  curves,  represented  by  the  vertical  bars,  are  the 
variance  added  due  to  the  turn.  The  slope  results  in  a  diffusion  coefficient  of  3.1210'10m2/s  and  the  jumps 
in  the  plot  correspond  to  turn-induced  variance  of  1.0110'8m2,  both  of  which  are  in  good  agreement  with 
Ramsey’s  experimental  results. 

Finally,  FEMLAB  can  be  used  to  verify  the  quadratic  variance  growth  seen  with  both  abrupt- 
transition  1/r  and  Laplacian  fields  by  using  the  geometry  shown  in  figure  26a  with  an  increasing  Peclet 
number.  The  results  in  figure  26b  are  for  Peclet  numbers  ranging  from  500  to  10,000,  where  the  diffusion 
coefficient  was  varied  while  the  geometry  was  held  constant. 
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Figure  26  -  (a)  Geometry  for  FEMLAB  results  in  figure  26b.  w=50pm, 

rc=  125  pm,  L^=  1000pm.  (b)  FEMLAB  variance  results  for  the  multi-turn 
channel  geometry  in  figure  26a  for  Pe=500  to  10,000. 

Figure  26b  shows  the  growth  rate  approaching  the  quadratic  infinite  Peclet  growth  as  the  Peclet 
number  is  increased.  It  also  shows  the  effects  of  the  other  regimes  on  the  growth  rate,  where  diffusion  is 
more  prevalent.  As  diffusion  is  increased,  it  shows  the  movement  from  the  high  Peclet  super-linear  to 
lower  Peclet  linear  and  sub-linear  growth.  While  the  high  Peclet  super-linear  growth  is  well  understood 
to  be  due  to  the  quadratic  dependence  on  the  number  of  turns  traversed,  as  seen  in  (a  19)  and  (a24)  in  the 
Appendix,  the  reason  for  the  trend  to  sub-linear  growth  is  not  as  apparent.  This  trend  in  the  growth  rate 
involves  the  interaction  of  field  and  geometry  based  dispersion,  axial  diffusion,  and  transverse  diffusion, 
which  all  interact  differently  as  the  Peclet  number  is  changed.  At  first,  this  transition  from  super-linear  to 
sub-linear  might  be  from  the  measurement  limitations  mentioned  in  section  3.4.  However  this  is  not  the 
case.  Such  measurement  limitations  come  from  the  fact  that  it  is  impossible  for  the  variance  to  grow 
beyond  the  square  of  the  length  of  the  straight  section  where  the  band  variance  is  measured.  However, 
in  this  case  the  late-time  asymptote  of  such  a  limitation  in  the  measurement  is  at  1-10'6  m2,  which  is 
several  orders  of  magnitude  above  the  current  set  of  results,  and  is  therefore  not  limiting  the  growth  of  the 
variance  measurements.  It  is  interesting  to  note,  however,  that  all  regimes  will  eventually  approach  the 
late-time  asymptote,  and  it  is  something  that  should  not  be  ignored  when  making  such  measurements. 
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3.7  High-Pure  Advection  Analytical  Dispersion  with  0-Assumption  Field  Model 
In  this  section,  some  of  the  results  gathered  from  the  analytical  abrupt-transition  1/r  field  variance 
and  particle  model  Laplacian  field  variance  simulations  in  the  previous  section  will  be  used  to  obtain  a 
high-Peclet  analytical  description  with  the  full  Laplacian  transition  field  structure. 

Since  a  first-principles  derivation  is  not  possible  to  develop  the  analytical  expression  for  the  full 
Laplacian  field  variance,  the  relationship  between  the  abrupt-transition  and  Laplacian  field  variances 
needs  to  be  characterized.  From  this  characterization  a  relational  factor  dependent  on  the  problem  space 
can  be  defined  to  correct  the  abrupt-transition  field  analytical  variance  descriptions  to  form  a  complete 
Laplacian  analytical  description  of  the  band  variance. 

From  section  3.5  it  is  known  that  there  is  a  constant  factor  between  the  two  descriptions  for  any 
number  of  turns  for  a  given  w/rc.  So,  the  only  thing  that  needs  to  be  characterized  is  the  factor  between 
the  two  expressions  for  a  w/rc  range  for  even-  and  odd-numbered  turns.  To  do  this,  a  series  of  particle 
simulations  were  run  to  determine  the  Laplacian  variance  results,  which  were  then  compared  with  the 
analytical  expression  for  the  abrupt-transition  field  variance.  The  factors  for  even-  and  odd-  numbered 
turns  are  shown  in  Figures  27a  and  27b. 
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Figure  27  -  Correction  factors  for  the  abrupt-transition  field  variance  equations  in  (30)  for  even-numbered  turns,  (a), 
and  odd-numbered  turns,  (b). 


Using  the  correction  factors  in  figure  27,  which  are  weak  functions  of  w/rc,  a  new  expression  for 
the  variance  using  the  Laplacian  transition  field  structure  can  be  created,  and  this  is  shown  here  with  an 
expanded  version  in  the  Appendix: 
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a  =  mod(«,2)(-0. 12945)  +  (1  -mod(n,2))(-0.23867) , 
P  =  mod(«, 2)(0.02408)  +  (1  -  mod(n,2))(-0.83142) , 
y  =  mod(n, 2)(0.99705)  +  (1  -  mod(n,2))(4.28054) 
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Laplacian  Particle  Simulations  Variance  Surface 


Laplacian  Analytical  Variance  Surface 
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Figure  28  -  A  surface  representing  the  variance  of  a  band  after  traveling  through  a  channel,  with 
w/rc  ranging  from  0.1  to  1.0  and  n  ranging  even  numbered  from  2  to  20,  as  described  by  the 
particle  simulation  in  (a)  and  (35)  in  (b). 
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Figure  29  -  The  error  between  the  two  surfaces  seen  in  figure  28. 
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To  verify  the  new  expression  in  (35),  a  surface  was  created  over  a  parameter  space  consisting  of 
w/rc  and  the  number  of  turns  using  the  particle  simulator.  The  number  of  turns  is  the  even-numbered 
turns  up  to  20  and  w/rc  ranges  from  0.1  to  1.0,  a  superset  of  the  typical  w/rc  design  space  for  serpentine 
channels  which  goes  up  to  0.4.  The  surfaces  are  seen  in  figure  28  and  the  error  between  the  two  is 
illustrated  in  figure  29. 

The  error  in  figure  29  is  computed  from  the  difference  of  the  two  surfaces  normalized  by  the 
surface  described  by  the  particle  simulation.  Over  nearly  the  entire  surface,  the  error  remains  less  than 
1%.  Only  in  the  corner  defined  by  w/rc  values  near  1.0  does  the  error  begin  to  increase  to  approximately 
6%.  The  source  of  this  slightly  increasing  inaccuracy  remains  to  be  investigated.  However,  since  this  is 
an  area  of  the  parameter  space  that  would  rarely,  if  ever,  be  used  for  useful  designs,  the  error  has  no  real 
effect  on  the  performance  of  the  model.  Further,  since  the  assumptions  appear  to  hold  for  the  entire 
feasible  design  space,  w/rc  equals  0.1  to  0.4,  the  performance  is  considered  acceptable.  If  needed,  the 
error  could  be  removed  by  renormalizing  the  analytical  variance  surface. 
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4.  Conclusion 


This  report  has  presented  a  new  model  of  the  high  Peclet  variance  from  the  electric  field  and 
channel  geometry  using  the  full  structure  of  the  electric  fields  described  by  complete  solutions  to 
Laplace’s  equation.  This  is  a  functional  form,  as  described  in  [6],  that  can  be  used  as  a  component  model 
in  a  separation  system  when  combined  with  other  subsystems  such  as  the  injector  and  detector.  As  a 
closed-form  algebraic  expression  this  model  can  be  used  for  rapid  computation  and  gradient-based 
optimization. 

To  create  this  model,  this  report  first  presented  additional  solutions  to  the  field  structure  in  the 
serpentine  turns,  which  gradually  relaxed  the  two  main  assumptions  of  linear  fields  in  the  curves  with  an 
abrupt  transition  to  the  straight  sections  made  in  [2, 3, 4, 5].  The  first  field  solution  presented  in  chapter  2, 
removed  the  assumption  of  linear  field  structure  using  a  1/r  field  in  the  curves  while  still  using  an  abrupt 
transition  to  the  straight  sections.  The  second  field  solution  presented  in  chapter  2  used  complete 
solutions  to  Laplace’s  equation  to  have  a  continuous  transition  between  sections  with  the  correct  structure 
in  straight  and  curved  sections,  as  verified  with  finite  element  software.  In  chapter  3,  an  analytical  model 
of  the  variance  of  a  band  traveling  through  the  serpentine  channels  was  created  using  the  abrupt-transition 
1/r  1 -assumption  field  structure.  Using  this  variance  model,  new  effects  were  predicted,  such  as 
incomplete  unskewing  after  complementary  turns  and  quadratic  growth  of  the  variance  after  n  even  or 
odd  turns,  which  were  confirmed  in  FEMLAB  simulations.  Finally,  this  abrupt-transition  model  was 
combined  with  results  from  the  Matlab  particle  simulations  using  the  O-assumption  full  Laplacian  field 
structure  to  create  the  complete  high  Peclet  variance  model. 

4.1  Future  Work 

In  its  current  functional  form,  the  high  Peclet  variance  model  will  be  combined  with  the  method 
of  moments,  M.O.M.,  model  developed  at  CMU  [to  be  published]  to  describe  the  variance  produced  by  a 
serpentine  of  any  number  of  turns  in  any  regime  of  figure  13.  Ultimately,  the  variance  descriptions  and 
insight  gained  from  the  high  Peclet  variance  expressions  in  this  work  will  be  used  to  create  a  component 
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based  model  for  serpentine  channels,  which  will  be  useful  for  creating  any  general  serpentine  geometry 
with  uniform  width.  This  component  model  will  be  only  for  the  turns  in  serpentine  channels  combined 
with  the  straight  section  model  of  the  which  will  allow  general  serpentines  to  be  constructed 

piece  by  piece.  The  challenge  in  converting  the  functional  variance  descriptions  in  this  report  will  be 
modifying  them  to  accept  various  initial  conditions,  since  the  input  to  the  component  model  will  be  the 
output  of  the  previous  straight  section  after  an  arbitrary  number  of  turns,  rather  than  a  well-defined 
injector. 
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Appendix 


A.l  Buckingham  7C-Theorem 

Background:  See  Harald  Hanche-Olsen’s  notes  in  [15]  for  further  details.  Provided  below  are  a 
partial  proof  and  examples  of  the  Buckingham  71-Theorem  from  previous  sections  of  this  report: 

Assuming  variables  with  values  and  units  that  describe  a  physical  function,  0(Rh...,Rn), 
according  to  the  following  notation: 

Rj=v(Rj  )[/?.],  (a5) 

where  Rj  is  the  variable  with  value  v(Rj)  and  units  [Rj\ 

Further,  it  is  assumed  the  units  of  the  variable(s)  can  be  broken  down  into  fundamental 
units  such  as  m,  kg,  s,  K,  using  a  consistent  set  such  as  MKS.  These  units  can  be  represented  as 
F1,...,Fm.  This  allows  representation  of  the  units  of  the  variable(s)  as  a  product  of  powers: 

m 

=  0 1  =  l.-,n),  (a6) 

i= 1 

where  it  is  necessary  that  the  units  themselves  are  independent  such  that: 

m 

IP  =!  (*,=..  =  *„  =0)  (a7) 

i= 1 

Now  a  matrix  of  the  powers  in  the  product  series  is  constructed,  whose  dimensions 
depend  on  the  number  of  variables,  n ,  and  the  number  of  fundamental  units  used,  m. 

an 

aml 

Using  these  properties  we  can  define  a  non-dimensional  group  from  the  original 
variables.  Such  a  group  is  constructed  out  of  a  series  of  products  of  the  original  variables  of  the 
system,  Rj,  raised  to  the  necessary  powers,  A n.  Such  a  group  would  have  units  defined  as: 

m 

[R*  ■  ■  ■  ]  =  Yl  ip+  "+fliA  (a9) 

i= 1 
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For  this  group  of  units  to  be  dimensionless,  or  for  (a9)  to  equal  one,  the  power  of  the 
product  series  must  equal  zero,  which  is  equivalent  to  AA^=  0,  where  >yT.  Therefore, 

finding  the  null-space  of  A,  N(A),  will  provide  the  set  of  powers  necessary  to  have  a  set  of  non- 
dimensional  variables.  It  can  provide  the  maximal  set  of  non-dimensional  variables  \nx ,. .  .,^_r(A)] 
without  knowledge  of  the  underlying  physical  equations  that  are  the  source  of  the  dimensional 
analysis.  r(A)  is  the  rank  of  A,  therefore,  n-r(A)  is  the  number  of  free  variables,  which  set  the 
number  of  independent  vectors  in  the  null  space  of  A  for  our  under-specified  linear  system.  The 
vectors  that  compose  N(A)  then  define  the  form  of  the  non-dimensional  variables. 


Example  1  (from  section  3.1) 


An  example,  based  on  the  non-dimensionalization  of  the  function  describing  the  variance 
of  a  band  traveling  through  serpentine  channels  in  any  of  the  four  regimes  shown  in  figure  13 
begins  with  the  function  describing  the  variance: 

<&(w,rc,lbt,Ec,D,ju)  (alO) 


where  w  is  the  channel  width,  rc  is  the  center  radius  of  the  turn,  Lbt  is  the  length  between  turns,  Ec 
is  the  constant  field  strength  in  the  straight  sections,  D  is  the  diffusion  coefficient,  and  //  is  the 
mobility  of  the  separation  species.  Each  of  these  variables,  Rj,  contributes  to  the  set  of 
fundamental  units  Fh  to  create  the  set  of  F=  { m,kg,s,A } .  This  allows  the  matrix  of  variable  unit 
powers  A  to  be  constructed: 


1111  2  0 

0  0  0  1  0  -1 

0  0  0  -3  -1  2 

0  0  0  -1  0  1 


(all) 


where  the  columns  represent  variables,  [w,rc,LbhEc,D,ju\,  and  the  rows  represent  fundamental 
units,  [m,kg,sA]T- 

Finding  the  nullspace  of  A,  A(A)  provides  the  following  result: 
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N(A)  = 


-1  -1  1 

1  0  0 

0  1  0 

0  0  1 

0  0-1 

0  0  1 


(a  12) 


The  powers  of  the  variables  in  each  of  the  non-dimensional  groups  are  represented  in  the 
chosen  basis  of  N(A).  For  this  example,  the  three  non-dimensional  groups  are  therefore: 

rc  h,  MEcw 


K  = 


W  W 


D 


(al3) 


Example  2  (from  Section  3.5) 

Following  the  same  procedure  for  the  analysis  of  high-pure  advection  regime  results  in 
the  following  trivial  example: 

00,rc),  (a  14) 

where  the  triviality  lies  in  the  fact  that  there  are  no  time-dependent  parameters,  since  diffusion  is 
not  an  issue.  Therefore  the  variance  fundamentally  only  depends  on  w  and  rc.  This  gives  the 
following  results  for  the  non-dimensional  parameters  of  which  there  will  be  one,  n-r(A)= 1: 

A  =  [l  l]  (al5) 


N(A)  = 


-1 

1 


(a  16) 


K  - 


w 


A.2  Expanded  Algebraic  Form  of  Equation  30. 

Expanding  Equation  30  from  the  text,  repeated  here: 


_2  _ 
&  ab  ~ 


k1  In2 (b / a)  V 


w 


2.*,,  \  2  (V+^2)((w+u)3-«3)^ 


Nxr  +  N 2{w  +  2a  —  r)  -■ 


a  \ 


3w 


(al7) 


dr 


(al8) 
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creates  the  following  closed  form  algebraic  for  the  high-pure  advection  variance  using  the  abrupt 
1/r  fields  in  the  turns: 


a 


2 

ab 


n 2  In  2(bl  a) 
15  w2 


'l5G2  +l5a4(Nl+N2)2  +30a3w(N1+N2)2  + 
l0a2w2(3N2  +4NlN2+3N22)  +  5aw3(3N2  +2NlN2+3N22)  + 
w\3N2  +NlN2+3N2)-10G{Nl+N2)(3a2  +3aw  +  w2) 


,(al9) 


where 


G  =  —(Nl+N2  XO  +  a)3  -  a3 ) 
3  w 


Nt  = 


n  +  mod(n,2) 


(a20) 


(a21) 


n2  = 


n  +  mod(n  + 1,2)  - 1 


(a22) 


A.3  Expanded  Algebraic  Form  of  Equation  35. 

The  expanded  form  of  35  is  identical  to  (a  19)  except  for  the  correction  factor  in  (36), 
which  results  in  the  following: 


2  7T2  \n2  (b  /  a) \(  2  ,  Ar  ,  ,2  (N1+N2)((w  +  ay -a3) 

°,ap  =  4 - 3 - J  Ns  +  N2(w  +  2a  -  r) - 


w 


2  .  K2  ln2(b/a) 

^lap  ~  A 


3  w 


r\5G2  +\5a4(Nl+N2)2  +30a3w(Nl  +N2 )2  + 


dr 


(a23) 


J 


15w 


10a2w2(3Aj  +4N1N2+3N2  )  +  5aw3(3N1  +2NlN2+3N2  )  + 


,(a24) 


(3Nl  +  NlN2  +  3N2  )-  10G(N1  +  N2){3a  +  3aw  +  w  ) 

where  G,  Nu  and  N2  are  the  same  as  in  (a20),  (a21),  and  (a22)  respectively,  and  from  (36),  (37), 
(38),  and  (39): 


(a25) 

(a26) 

(a27) 

(a28) 


(  \ 

2 

(  \ 

A  =  a 

w 

+  P 

W 

vrc; 

Vrcy 

+  r , 

a  =  mod(n,2)(-0. 12945)  +  (1  -  mod(n,2))(-0.23867) , 
ft  =  mod(n,2)(0.02408)  +  (1  -  mod(n,2))(-0.83142) , 
y  =  mod(n,2)(0.99705)  +  (1  -  mod(n,2))(4.28054) 
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Symbol  Glossary 


a  inner  turn  radius 

b  outer  turn  radius 

c  concentration 

cavg  average  concentration 

Cn  series  coefficients  for  the  Laplacian  solutions  in  the  curved  sections 

D  diffusion 

Dn  series  coefficients  for  the  Laplacian  solutions  in  the  straight  sections 

Excurve  transverse  electric  field  component  in  the  straight  section  (Laplacian  transition  field) 

Ey  axial  electric  field  component  in  the  straight  section  (abrupt-transition  1/r  field) 

E'y straight  axjai  electric  field  component  in  the  straight  section  (Laplacian  transition  field) 

Ercurve  radial  electric  field  component  in  the  curved  section  (Laplacian  transition  field) 

Ec/)  circumferential  electric  field  component  in  the  curved  section  (abrupt-transition  1/r  field) 

E^urve  circumferential  electric  field  component  in  the  curved  section  (Laplacian  transition  field) 

G  common  factor  found  in  the  abrupt  and  Laplacian  analytical  variance  equations 

k  Boltzmann’s  constant 

k  wave  number  in  the  Laplacian  series  solutions 

/  straight  channel  length  for  uniform  symmetric  single  turns 

Lbt  length  between  turns  for  serpentine  channels 

Ldet  length  of  the  last  leg  of  serpentine  going  to  the  detector 

Linj  length  of  first  leg  of  serpentine  coming  from  the  injector 

n  number  of  turns  in  the  serpentine  channel 

n±  ion  concentration  in  the  double  layer  for  positive(+)  and  negative(-)  ions 

n±°  ion  concentration  in  the  bulk  for  positive(+)  or  negative(-)  ions 

Nj  indexing  modulo  function  for  odd  turns  found  in  analytical  variance  expressions 

N2  indexing  modulo  function  for  even  turns  found  in  analytical  variance  expressions 

Pe  Peclet  number  based  on  width 

pts  number  of  particles  in  the  Matlab  particle  simulation 

qe  charge  of  an  electron 

r  radial  position  in  the  curved  section  of  the  channels 

rc  mean  channel  width 

skew  the  skew  of  a  band  traveling  in  a  linear  electric  field  in  a  turn 
t  time 

T  temperature 

Tbt  travel  time  ratio  of  width  to  length  between  turns 

Ti  time  for  a  particle  to  travel  the  inner  most  radius  in  the  turn 

T0  time  for  a  particle  to  travel  the  outer  most  radius  in  the  turn 

Tt  travel  time  ratio  of  width  to  length  of  turn  along  center  radius 

u  channel  velocity  field 

u0  mean  channel  velocity 

V0  applied  electric  potential 

w  channel  width 

v  transverse  direction  in  the  straight  section  of  the  channels 

y  axial  direction  in  the  straight  section  of  the  channels 

z  ion  valency 

a  component  of  the  correction  factor  X 

P  component  of  the  correction  factor  X 

£  buffer  permittivity 

<E>  electric  potential 
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(f)  angular  position  in  the  curved  section  of  the  channels 

y  component  of  the  correction  factor  X 

T)  Laplacian  series,  coefficient  coupling  block  matrix  component 

X  correction  factor  for  the  abrupt  analytical  variance  to  produce  the  Laplacian  analytical  variance 

//  species  mobility 

jUv  viscosity  of  the  buffer 

0  turn  angle 

p  free-charge  density 

(jab  abrupt-transition  Hr  field  analytical  variance 

< 7iaP 2  full  Laplacian  transition  field  analytical  variance 

T  Laplacian  series,  coefficient  coupling  block  matrix  component 

£  Laplacian  series,  coefficient  coupling  block  matrix  component 

£  zeta  potential 


Normalized  Variables 

c  concentration,  c-dc0 

t  time,  t  =t/(w/u0 ) 

u  velocity,  u  -u!u0 

v*  transverse  position,  x*=x/w 

y*  axial  position,  y*=y/w 
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